/*
** Xin YUAN, 2024, BSD (2)
*/

////////////////////////////////////////////////////////////////////////////////

//triangle

// bary center method
// P - V1 = (V2 - V1) * u + (V3 - V1) * v
template <typename T>
inline bool Triangle_IsPointIn(const Vec3<T>& p, const Vec3<T>& v1, const Vec3<T>& v2, const Vec3<T>& v3,
				T& u, T& v) noexcept
{
	Vec3<T> const ep(p - v1);
	Vec3<T> const e2(v2 - v1);
	Vec3<T> const e3(v3 - v1);
	T d = e2.coord.x * e3.coord.y - e2.coord.y * e3.coord.x;
	if ( DefaultCompareTrait<T>::IsEQ(d, T(0)) )
		return false;
	d = T(1) / d;
	u = (ep.coord.x * e3.coord.y - ep.coord.y * e3.coord.x) * d;
	v = (ep.coord.y * e2.coord.x - ep.coord.x * e2.coord.y) * d;
	return DefaultCompareTrait<T>::IsGE(u, T(0))
		&& DefaultCompareTrait<T>::IsGE(v, T(0))
		&& DefaultCompareTrait<T>::IsLE(u + v, T(1))
		&& DefaultCompareTrait<T>::IsEQ(u * e2.coord.z + v * e3.coord.z, ep.coord.z)
		;
}

//plane

// N dot P + d = 0
// plane : (Nx, Ny, Nz, d)

template <typename T>
inline T Plane_Dot(const Vec4<T>& p, const Vec4<T>& v) noexcept
{
	return p.Dot(v);
}
template <typename T>
inline T Plane_DotCoord(const Vec4<T>& p, const Vec3<T>& v) noexcept
{
	return p.plane.a * v.coord.x + p.plane.b * v.coord.y + p.plane.c * v.coord.z + p.plane.d;
}
template <typename T>
inline T Plane_DotNormal(const Vec4<T>& p, const Vec3<T>& v) noexcept
{
	return p.plane.a * v.coord.x + p.plane.b * v.coord.y + p.plane.c * v.coord.z;
}

template <typename T>
inline Vec4<T> Plane_FromPointNormal(const Vec3<T>& p, const Vec3<T>& n) noexcept
{
	return Vec4<T>(p.coord.x, p.coord.y, p.coord.z, -p.Dot(n));
}
template <typename T>
inline Vec4<T> Plane_FromPoints(const Vec3<T>& v1, const Vec3<T>& v2, const Vec3<T>& v3) noexcept
{
	return Plane_FromPointNormal(v1, (v2 - v1).Cross(v3 - v1));
}

template <typename T>
inline Vec3<T> Plane_GetNormal(const Vec4<T>& p) noexcept
{
	return Vec3<T>(p.plane.a, p.plane.b, p.plane.c);
}

template <typename T>
inline Vec4<T> Plane_Normalize(const Vec4<T>& p) noexcept
{
	T const n = Plane_GetNormal(p).Length();
	return p / n;
}

//reflection (mirror)

template <typename T>
inline Mat4x4<T> Mat_Reflection(const Vec4<T>& plane) noexcept
{
	return Mat4x4<T>(
			//col 0
			-T(2) * plane.plane.a * plane.plane.a + T(1),
			-T(2) * plane.plane.a * plane.plane.b,
			-T(2) * plane.plane.a * plane.plane.c,
			-T(2) * plane.plane.a * plane.plane.d,
			//col 1
			-T(2) * plane.plane.b * plane.plane.a,
			-T(2) * plane.plane.b * plane.plane.b + T(1),
			-T(2) * plane.plane.b * plane.plane.c,
			-T(2) * plane.plane.b * plane.plane.d,
			//col 2
			-T(2) * plane.plane.c * plane.plane.a,
			-T(2) * plane.plane.c * plane.plane.b,
			-T(2) * plane.plane.c * plane.plane.c + T(1),
			-T(2) * plane.plane.c * plane.plane.d,
			//col 3
			T(0), T(0), T(0), T(1));
}

//shadow

// light.coord.w : 0 --- directional light, 1 --- point light
template <typename T>
inline Mat4x4<T> Mat_Shadow(const Vec4<T>& light, const Vec4<T>& plane) noexcept
{
	T const d = -plane.Dot(light);
	return Mat4x4<T>(
			//col0
			plane.plane.a * light.coord.x + d,
			plane.plane.b * light.coord.x,
			plane.plane.c * light.coord.x,
			plane.plane.d * light.coord.x,
			//col1
			plane.plane.a * light.coord.y,
			plane.plane.b * light.coord.y + d,
			plane.plane.c * light.coord.y,
			plane.plane.d * light.coord.y,
			//col2
			plane.plane.a * light.coord.z,
			plane.plane.b * light.coord.z,
			plane.plane.c * light.coord.z + d,
			plane.plane.d * light.coord.z,
			//col3
			plane.plane.a * light.coord.w,
			plane.plane.b * light.coord.w,
			plane.plane.c * light.coord.w,
			plane.plane.d * light.coord.w + d);
}

//ray

// P = P0 + tU

template <typename T>
inline bool Ray_IntersectPlane(const Vec3<T>& s, const Vec3<T>& u, const Vec4<T>& p, T& t) noexcept
{
	Vec3<T> const n(Plane_GetNormal(p));
	T const d = n.Dot(u);
	//if direction and plane-normal are not perpendicular
	if ( DefaultCompareTrait<T>::IsEQ(d, T(0)) )
		return false;
	t = (-p.plane.d - n.Dot(s)) / d;
	return true;
}

// bounding sphere : vec4<T> : (x0, y0, z0, radius)

template <typename T>
inline bool Ray_IntersectSphere(const Vec3<T>& s, const Vec3<T>& u, const Vec4<T>& p,
				T& t1, T& t2) noexcept
{
	Vec3<T> const v(s - Vec3<T>(p.sphere.x, p.sphere.y, p.sphere.z));
	T A = u.Dot(u);
	T const B = T(2) * u.Dot(v);
	T D = SQUARE(B) - 4 * A * (v.Dot(v) - SQUARE(p.sphere.r));
	if ( DefaultCompareTrait<T>::IsLT(D, T(0)) )
		return false;
	D = FloatTraits<T>::sqrt(D);
	A = T(1) / (T(2) * A);
	t1 = (-B + D) * A;
	t2 = (-B - D) * A;
	return true;
}

//pick

// to view coordinates

template <typename T>
inline void Pick_FromViewport(const T& x, const T& y,
				const Mat4x4<T>& vp, const Mat4x4<T>& proj,
				Vec2<T>& v) noexcept
{
	T const w_2 = vp.v._11;
	T const h_2 = -vp.v._22;
	v.coord.x = ((x - (vp.v._41 - w_2)) / w_2 - T(1)) / proj.v._11;
	v.coord.y = (-(y - (vp.v._42 - h_2)) / h_2 + T(1)) / proj.v._22;
}

template <typename T>
inline void PickRay_OrthoRH(const Vec2<T>& p, const T& zn, Vec3<T>& s, Vec3<T>& u) noexcept
{
	s.Set(p.coord.x, p.coord.y, -zn);
	u.Set(T(0), T(0), T(1));
}
template <typename T>
inline void PickRay_OrthoLH(const Vec2<T>& p, const T& zn, Vec3<T>& s, Vec3<T>& u) noexcept
{
	s.Set(p.coord.x, p.coord.y, zn);
	u.Set(T(0), T(0), T(1));
}

template <typename T>
inline void PickRay_PerspectiveRH(const Vec2<T>& p, const T& zn, Vec3<T>& s, Vec3<T>& u) noexcept
{
	s.Set(T(0), T(0), T(0));
	u.Set(p.coord.x, p.coord.y, -zn);
}
template <typename T>
inline void PickRay_PerspectiveLH(const Vec2<T>& p, const T& zn, Vec3<T>& s, Vec3<T>& u) noexcept
{
	s.Set(T(0), T(0), T(0));
	u.Set(p.coord.x, p.coord.y, zn);
}

////////////////////////////////////////////////////////////////////////////////
